; mv idl.ps ~/Dropbox/Tex/paper.nssl/fig/theprof.eps
if !d.name eq 'PS' then begin
    device,xsize=22,ysize=12,yoffset=3
    !p.thick=2 & !x.thick=1 & !y.thick=1
end
!x.margin=[4,3] & !y.margin=[4,3]
!p.multi=[0,2,1]

dir1='/home/guerrero/eumag/nt-28-rf-movt/'
dir2='/home/guerrero/eumag/tac-28-rf-movt/'

file1a=dir1+'velocity.sav'
restore,file1a
print,'itac = 0'
itac=0
@rhoprof
print,rh
m=(size(rurms2))[2]

ithp=where(th ge -70 and th le -50)
nthp=(size(ithp))[1]
omp=total(omega[*,ithp],2)/double(nthp)
dlnompdr=deriv(r,alog(omp))

ithe=where(th ge -30 and th le -10)
nthe=(size(ithe))[1]
ome=total(omega[*,ithe],2)/double(nthe)
dlnomedr=deriv(r,alog(ome))

plot,r[1:45],dlnompdr[1:45],xrange=[0.61,0.96],charsize=1.5,xstyle=1, $
	yr=[-1.5,1.2],ystyle=1,thick=5,xtitle='!8r/R!6'+sunsymbol(),$
	title='a) !9d!D!8r!N!6ln!4X!6',/nodata,xticks=4
oplot,r[1:45],smooth(dlnompdr[1:45],2),li=2,thick=10,color=cgColor('red')
oplot,r[1:45],smooth(dlnomedr[1:45],2),li=0,thick=10,color=cgColor('red')

file2a=dir2+'velocity.sav'
restore,file2a
print,'itac = 1'
itac=1
@rhoprof
m=(size(rurms2))[2]

ithp=where(th ge -70 and th le -50)
nthp=(size(ithp))[1]
omp=total(omega[*,ithp],2)/double(nthp)
dlnompdr=deriv(r,alog(omp))

ithe=where(th ge -30 and th le -10)
nthe=(size(ithe))[1]
ome=total(omega[*,ithe],2)/double(nthe)
dlnomedr=deriv(r,alog(ome))

oplot,r[1:62],smooth(dlnompdr[1:62],2),li=2,thick=10,color=cgColor('blue')
oplot,r[1:62],smooth(dlnomedr[1:62],2),li=0,thick=10,color=cgColor('blue')
oplot,[0.6,0.98],[0.,0.],li=1,thick=2
oplot,[0.68,0.68],[-40,40],li=1,thick=2
oplot,[0.77,0.77],[-40,40],li=1,thick=2

oplot,[0.7,0.73],[-1,-1],li=2,thick=10,color=cgColor('red')
oplot,[0.73,0.76],[-1,-1],li=2,thick=10,color=cgColor('blue')
xyouts,0.78,-1.02,'!6pole',charsize=1.5
oplot,[0.7,0.73],[-1.3,-1.3],li=0,thick=10,color=cgColor('red')
oplot,[0.73,0.76],[-1.3,-1.3],li=0,thick=10,color=cgColor('blue')
xyouts,0.78,-1.302,'!6equator',charsize=1.5
;;;;;;;;;;;;;;
;;;;;;;;;;;;;;
;;;;;;;;;;;;;;

file1a=dir1+'velocity.sav'
restore,file1a
print,'itac = 0'
itac=0
@rhoprof
m=(size(rurms2))[2]
thg=where(th GE -84 AND th LT 84)

ibase=where(r ge 0.72 and r le 77)
rbmean=mean(r[ibase])
nbase=(size(ibase))[1]
omb=total(omega[ibase,*],1)/double(nbase)
dombdth=deriv(th,alog(omb))/rbmean

itop=where(r ge 0.90 and r le 95)
rtmean=mean(r[itop])
ntop=(size(itop))[1]
omt=total(omega[itop,*],1)/double(ntop)
domtdth=deriv(th,alog(omt));/rtmean

plot,th[thg],dombdth[thg],xrange=[-90,90],charsize=1.5,xstyle=1, $
	yr=[-0.02,0.02],ystyle=1,thick=5,xtitle='!890-!7h!6',$
	title='b) !8r!U-1!N!9d!D!7h!N!6ln!4X!6',/nodata,xticks=4
oplot,th[thg],dombdth[thg],li=2,thick=10,color=cgColor('red')
oplot,th[thg],domtdth[thg],li=0,thick=10,color=cgColor('red')

file2a=dir2+'velocity.sav'
restore,file2a
m=(size(rurms2))[2]
;
print,'itac = 1'
itac=1
@rhoprof
ibase=where(r ge 0.72 and r le 77)
rbmean=mean(r[ibase])
nbase=(size(ibase))[1]
omb=total(omega[ibase,*],1)/double(nbase)
dombdth=deriv(th,alog(omb))/rbmean

itop=where(r ge 0.90 and r le 95)
rtmean=mean(r[itop])
ntop=(size(itop))[1]
omt=total(omega[itop,*],1)/double(ntop)
domtdth=deriv(th,alog(omt))/rtmean

oplot,th[thg],dombdth[thg],li=2,thick=10,color=cgColor('blue')
oplot,th[thg],domtdth[thg],li=0,thick=10,color=cgColor('blue')
;
oplot,[-70,-45],[-0.01,-0.01],li=2,thick=10,color=cgColor('red')
oplot,[-45,-20],[-0.01,-0.01],li=2,thick=10,color=cgColor('blue')
xyouts,-15,-0.0105,'!6bottom',charsize=1.5
oplot,[-70,-45],[-0.015,-0.015],li=0,thick=10,color=cgColor('red')
oplot,[-45,-20],[-0.015,-0.015],li=0,thick=10,color=cgColor('blue')
xyouts,-15,-0.0155,'!6top',charsize=1.5

oplot,[-90,90],[0,0],li=1,thick=2
oplot,[0,0],[-1,1],li=1,thick=2

;;;;;;
;;;;;;
;;;;;;


;oplot,[0.85,0.9],[340,340],li=0,thick=10,color=cgColor('red')
;xyouts,0.91,300,'CZ02',charsize=2
;oplot,[0.85,0.9],[120,120],li=2,thick=10,color=cgColor('blue')
;xyouts,0.91,100,'RC02',charsize=2

end
